Gamma Calibration

Introduction and loading data

Detector response to incident gamma rays across the detector varies as a result of various effects at the channel and cube level, including the cube light yield, the channel coupling to the SiPMs, and attenuation effects within the WLS fibres. Due to the directional dependence of effective dose, it is key to ensure that the detector response is isotropic. As a result, we require calibration of these channel effects such that we measure the same distribution of ADC counts per channel for the same energy of incident particles, as otherwise rotation of the detector would change the measured effective dose for the same radiation field.

In 2020, data was taken using a Co\(^{60}\) source in a lab setting at Imperial. This was used to perform a channel-to-channel light response calibration, using the Compton edge from the first photopeak of Co\(^{60}\) as a reference point.

The majority of the code used for this analysis can be found in GammaCalib.py, but key sections will be reproduced here.

import GammaCalib as GammaCalib

import numpy as np
import pandas as pd

%matplotlib inline

The analysis code is contained within the CouplingAnalysis class, which is instantiated with a timing window for determining coincidence between signals in the detector.

analysis = GammaCalib.CouplingAnalysis(100)

Individual data files can be loaded into the CouplingAnalysis object via the add_data method, which takes the ROOT file of the data run. The peak value over pedestal for each signal is calculated on a channel-by-channel basis by subtracting the channel baseline from the peak of the waveform. Coincidence calculation is then performed to prepare the datafile for analysis. Background data is also loaded for the purposes of background subtraction during the analysis.

data_file = '/home/nr1315/Documents/Project/GammaCalibration/Data/lab_2020-09-16_Co60_lowToT.root'
bkg_file = '/home/nr1315/Documents/Project/GammaCalibration/Data/lab_2020-09-16_bkg_lowToT.root'

analysis.add_data('/home/nr1315/Documents/Project/GammaCalibration/Data/lab_2020-09-16_Co60_lowToT.root')
analysis.add_data('/home/nr1315/Documents/Project/GammaCalibration/Data/lab_2020-09-16_bkg_lowToT.root')

Channel fitting

Now that the data has been loaded, it is instructive to visualise the distribution of peak values over pedestal for a single channel.

This shows the Compton energy distribution from the Co\(^{60}\) source incident in the detector, which comprises of Compton scattering from the two photopeaks of the source. However, peak broadening effects in the detector distort these distributions, resulting in the shape seen in this plot. For this analysis, we focus on the distribution from the first photopeak of Co\(^{60}\) at 1173.2 keV, with a Compton edge energy of 963.3 keV. The broadening effects give rise to a non-trivial distribution, but here it is approximated by a half-Gaussian function from the peak of the distribution onwards. The ADC value of this Compton edge can be extracted from the value of ADC at the 1\(\sigma\) width of the distribution. By fitting a half-Gaussian distribution to each channel, this value can be extracted from each channel.

The gen_df function of the CouplingAnalysis class performs this fit for each channel and calculates the chi-squared value of each fit, and returns this information in a Pandas DataFrame. If a background dataset is provided, background subtraction is done at this stage as well.

analysis.gen_df(dfile = data_file, binsize = 128, maxval = 6144, bkg = bkg_file, ToTcut = 20, corrected = False)

/home/nr1315/miniconda3/lib/python3.7/site-packages/numpy/core/fromnumeric.py:3373: RuntimeWarning: Mean of empty slice.
  out=out, **kwargs)
/home/nr1315/miniconda3/lib/python3.7/site-packages/numpy/core/_methods.py:170: RuntimeWarning: invalid value encountered in double_scalars
  ret = ret.dtype.type(ret / rcount)

This fit can be visualised using the individual_gauss_fit function, which produces a plot of the individual channel distribution alongside the half-Gaussian fit and the residuals of the fit. This function can be used to scroll through each of the channels. An example channel is shown below.

analysis.individual_gauss_fit(data_file)

After the initial round of fitting, the reduced \(\chi^2\) values are checked against a user-defined threshold to decide on the quality of the fit. Here the threshold is set as 2, so any fits with a \(\chi^2\) value greater than this are identifed and labelled as poor fits. Then, the mean fit parameters \(\mu_i\) and \(\sigma_i\) are calculated from the good fits, and these are used as initial parameter guesses for re-fitting the channels with poor initial fits. This is done using the refit_mean function, and an example of the effect this has can be seen below.

Following this, the correction factors are calculated for each channel. The correction factor for the channel \(i\) on FPGA \(a\) is defined as

\[f_i = \frac{E_i}{{\langle E_i \rangle}^a}, \qquad E_i = \mu_i + \sigma_i,\]

where \({\langle E_i \rangle}^a\) denotes the mean over channels across FPGA \(a\). The peak value over pedestal distribution for each channel \(v_i\) is then scaled using these correction factors, such that the corrected values \(v^\prime_i\) are define as

\[v^\prime_i = \frac{v_i}{f_i},\]

for channel \(i\). The corrected values are saved in the DataFrame under the “corrected” key.

After this, the gen_df function is used to recalculate the fit for the corrected data to check that the correction factors have scaled the distribution in the correct direction. In order to verify the effects of this scaling effect, we look at the distribution of \(E_i\) for both FPGAs, before and after applying the correction factors to the peak value distributions.

Following application of the correction factors, there is a much better agreement between the different channel values of \(E_i\) for channels across each FPGA.